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ABSTRACT 


Recent developments at several levels of 
statistical turbulence modeling applicable to 
aerodynamics are briefly surveyed. Emphasis is on 
examples of model improvements for transonic, 
two-dimensional flows. Experience with the devel- 
opment of these improved models is cited to sug- 
gest methods of accelerating the modeling process 
necessary to keep abreast of the rapid movement of 
computational fluid dynamics into the computation 
of complex three-dimensional flows. 

NOMENCLATURE 


x = space location vector 

y = coordinate normal to surface 

z = coordinate in transverse direction 

6* = modeling coefficient 

6* = displacement thickness 

= Kronecker delta 

e = dissipation rate of kinetic energy 

= dissipation rate of Reynolds stress, u | u j 
v = kinematic viscosity of fluid 


C el = modeling coefficient in production 
of dissipation rate 

C £ 2 - modeling coefficient in destruction 

term of dissipation rate 

C v = modeling coefficient in eddy 
viscosity formulation 

D k = diffusion of kinetic energy 

D^j = diffusion of Reynolds stress, u [ u j 

k = mean kinetic energy of turbulence 

K = von Karman constant, 0.4 

L = characteristic length scale 


Q ij 


= turbulent eddy viscosity 

= fluid density 

= wall shear stress 

= pressure/rate-of-str ain production 
of Reynolds stress, u!u^ 

= specific dissipation rate, e/e*k 

= vorticity, 


a = (2a a ,i/z 

2 \ 9x j 3x i/ mn mn 


() = mean or averaged quantity 


= production of kinetic energy 


()' = turbulent fluctuating quantity 



t 


= production of Reynolds stress, 


= strain rate, 


3U . 
1 

3X j 



= time 


U i U j 


() i = vector component in ith direction, 

i = 1, 2, or 3 

()|j = ijth tensor 


u = velocity component in streamwise 
direction 


INTRODUCTION 


U e = velocity at boundary-layer edge 

v = velocity component normal to 
surface 

V = characteristic velocity scale 

w = velocity component in transverse 

direction 

x = coordinate in streamwise direction 


Numerical codes for the computation of fluid 
dynamical problems contain several distinct ele- 
ments. These include the basic equation set, 
grid-generation and nesting techniques for complex 
configurations, efficient solution algorithms, and 
appropriate boundary condition treatments. In 
addition, if the fluid flow is turbulent, the 
effects of all the scales of the flow field that 
cannot be resolved by the numerical scheme either 
in time or space must be modeled. It is through 
this "turbulence modeling" that the numerically 
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unresolved physics of fluid motion is introduced 
into the flow-field computation. 

The purpose of this paper, then, is to illus- 
trate the current status of "turbulence modeling" 
as it is applied to aerodynamic flow fields, and 
to suggest some means of accelerating progress in 
its development. That the latter is particularly 
important was recognized by two National Academy 
of Sciences (NAS) groups that studied the rate of 
progress of computational fluid dynamics (CFD) in 
general and turbulence modeling in particular. In 
1983, the committee headed by Robert Smelt con- 
cluded, "If the gap between progress in turbulence 
modeling and developments in numerical algorithms 
and computer hardware is allowed to continue, it 
could seriously impair the accuracy and usefulness 
of CFD towards the end of this century" (Smelt, 
1983). This view was reflected again in a yet to 
be published 1986 National Research Council 
Report. As turbulence is a feature of almost 
every aerodynamically important flow field, this 
situation deserves the serious attention of the 
CFD community. 

To examine why it is generally perceived that 
progress in turbulence modeling has been rather 
slow, it is important to understand what is 
involved in the creation or modification of a 
turbulence model. Basic to this understanding is 
the knowledge of the actual features of fluid 
turbulence and what the models can or cannot be 
expected to emulate. 

What then is fluid turbulence? Turbulence is 
a fluid flow, not a fluid property. And it pos- 
sesses many extremely complex characteristics. At 
a point in space, even when the gross flow is 
steady, the turbulent flow appears to be a chaoti- 
cally, or randomly, unsteady event about some mean 
value. Over an extended space, however, some 
coherence can be observed in the large-scale tur- 
bulent motions under certain conditions. The 
turbulent flow is always three-dimensional, even 
when it exists within a mean motion that is one- 
or two-dimensional. This three-dimensionality is 
very important, for it provides the degrees of 
freedom to allow local vortices to stretch, and in 
so doing it transmits energy from larger-scale 
motions to those that are smaller. As a conse- 
quence, the turbulent flow possess an extremely 
large range of scales in time and space, and the 
range of scales between the largest, related to 
body dimension, and the smallest, where the turbu- 
lence dissipates quickly, increases with larger 
Reynolds numbers. 

For an aircraft, the largest scales are a few 
percent of its length and carry most of the 
kinetic energy of the turbulence. The smallest 
scales, where the turbulence dissipates into heat, 
are about one thousandth of the local boundary- 
layer thickness, depending somewhat on the 
Reynolds number of the flow in question. Turbu- 
lence is also very diffusive and promotes mixing 
of sensible properties, such as temperature, 


chemical constituents, or momentum. For example, 
within aircraft engines this is a favorable prop- 
erty which enhances the distribution of fuel prior 
to combustion. On a surface of an aircraft, how- 
ever, this enhanced mixing contributes to 
increased drag. Finally, turbulence is dissipa- 
tive in that the fluid viscosity can destroy the 
motions of the smallest scales and convert their 
energy into heat. It has been observed experi- 
mentally that the large scales of turbulence are 
unique to individual flow situations, whereas the 
smaller scales tend to have a universal character. 
The range of scales broadens with increased 
Reynolds number. 

Up to about 15 years ago, the limitations of 
available computational power forced all turbu- 
lence calculations to be confined to the use of 
steady-state, statistical descriptions of the 
turbulence. Most of the features described in the 
previous paragraph could not be considered in 
detail. Turbulence models were applied only to 
averaged equations, and the information contained 
in the real dynamical features could only be 
restored, approximately, through the use of empir- 
ical constants that were contained in the statis- 
tical model formulations. Then, about 15 years 
ago, computers were developed that had sufficient 
power to allow the computation of the dynamic 
behavior of some very simple flow fields at very 
low Reynolds numbers where the range of scales 
present is relatively small. Underlying these 
computations is the assumption that the Navier- 
Stokes equations, together with the energy- 
transport equation, constitute the basic equations 
for continuum fluid flow. It is further assumed, 
without rigorous mathematical proof, that these 
three-dimensional nonlinear equations can be 
solved accurately in time and space and that the 
solutions lie within narrow bounds despite large 
numbers of uncertainties that exist in defining 
initial and time-dependent boundary conditions. 

Although these computations can capture all 
of the features of real turbulence described ear- 
lier, they are extremely costly and have been 
confined to studies of the "physics" of turbulence 
in idealized flow situations. These studies have 
resulted in descriptions of the dynamic mechanisms 
that take place in free and wall-bounded turbulent 
flows in much greater detail than can be measured 
in physical experiments. They also provide infor- 
mation on quantities that can only be measured 
inaccurately, if at all, for example, static pres- 
sure fluctuations or local rates of turbulence 
dissipation. Even though these computations deal 
with relatively simple flow fields of limited 
extent, they fill the capacity of the largest 
computers and require many tens of hours of com- 
puter time to reach statistically meaningful 
results. Some of these computations have yielded 
results that will be of inestimable value to the 
turbulence modeler, but this approach cannot be 
considered an engineering design tool even in the 
reasonably near future. The reasons for this are 
illustrated in figure 1. 
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Figure 1 compares the computer speed and 
memory requirements for different classes of aero- 
dynamic computational methods with the capabil- 
ities of some large existing computers and of some 
that are projected for the near term (Chapman, 
1978). For example, application of turbulence 
simulation where all the significant scales of 
turbulence are resolved for an airfoil (A, fig. 1) 
would require computer capabilities about 4 orders 
of magnitude larger than that labeled "Large Eddy 
Simulation." Clearly, such calculations cannot be 
considered in the foreseeable future. The regions 
labeled Large Eddy Simulation were estimated by 
Chapman under the assumption that only the sub- 
layer regions of boundary layers need to be 
resolved to the small scales and that the the 
outer parts of the boundary layers, the regions of 
separation, and shedding vortices could be treated 
inviscidly (Chapman, 1978). 

Even with the major economies afforded by 
these assumptions, the abilities to perform time- 
accurate simulations of the turbulent flow fields 
about aerodynamic bodies lie outside the abilities 
of the computers projected in the immediate 
future. Chapman’s estimates relative to the fully 
resolved computational needs, however, show the 
importance of the development of good sub-grid 
turbulence models to account for turbulence scales 
smaller than the grids employed for the bulk of 
the flow. In addition, if a good sub-grid model 
were developed for the near-wall region of the 
boundary layer, the computer requirements could be 
reduced further from those shown, and large eddy 
simulations of turbulent flow about an airfoil 
could be anticipated in the reasonably near 
future, not necessarily as a design technique, but 
at least as a research tool. The development of 
good sub-grid turbulence models for large-eddy 
simulations, then, is a research topic that should 
be given considerably more emphasis than it has 
had in the United States. The French and the 
Japanese have recognized this for several years, 
and are well ahead of us in this research topic. 

Figure 1 also illustrates that projections of 
computer requirements for the solution of turbu- 
lent aerodynamic flow fields about aircraft shapes 
with the Reynolds-averaged Navier-Stokes equations 
are within the capabilities of near future com- 
puters. In particular, the current NAS computer 
is capable of handling the turbulent flow over a 
complete wing with the Navier-Stokes approach and 
use of simple algebraic models of turbulence. It 
is the turbulence modeling for this class of com- 
putation that is the focus of the remainder of 
this paper. 


STATISTICAL TURBULENCE MODELING 


Reynolds-Averaged Equations 

Statistical turbulence modeling begins with 
the derivation of the equations to represent the 


mean or average motions of the turbulent flow. 
This process can be illustrated simply by examin- 
ing only the single equation representing the 
velocity parallel to the surface in an incompres- 
sible fluid, the instantaneous u-component equa- 
tion of the Navier-Stokes system of equations: 
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For a steady-state turbulent flow, the dependent 
variables are then expanded as the sum of their 
mean parts plus their turbulent fluctuations, 
namely, 


u ( x , t ) = u( x ) + u'(x,t) (2) 

where x represents a point in space and t is 
time. The overbar indicates an average quantity, 
and the prime represents the instantaneous fluctu- 
ation. The components of velocity v, normal to 
the surface, and w, in the transverse direction, 
are expanded similarly. Note that the steady- 
state requirement of this example eliminates time 
as an independent variable in the mean quantity. 


When these velocities are introduced into 
equation (1), and the resulting equation is aver- 
aged over time, there results the following: 
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It is noted that only the nonlinear moments of the 
fluctuating turbulent velocity components have 
survived the averaging process, which has elimi- 
nated most of the real physical details of the 
turbulence that were described earlier. The cha- 
otic time dependence, the phase relationships 
between the turbulence velocity components, and 
the multitude of turbulence scales are merged into 
just three averaged moments of the fluctuating 
velocity components. The resulting equation 
resembles the original Navier-Stokes equation, 
restricted to steady flow, except for the averaged 
moments of velocity fluctuations that reside in 
positions corresponding to those of the viscous 
stress terms and which, consequently, are called 
Reynolds stresses, after their original discoverer 
(Reynolds, 1895) . 

Note that these Reynolds stresses constitute 
three new dependent variables, whereas no new 
additional equations resulted from the averaging 
process. Thus, more unknowns are created than 
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equations, and the problem is not closed, at this 
level, unless the Reynolds stresses can be 
expressed in terms of the independent variables 
(which are the coordinates) or the dependent vari- 
ables (which are the velocity components) or both. 
If the six equations for the Reynolds stresses are 
derived from moments of the Navier-Stokes equa- 
tions, it is found that these new equations con- 
tain many more than six new higher-order terms, 
including those that contain averages of the fluc- 
tuating velocity derivatives (Chou, 19^5). The 
proliferation of dependent variables greatly 
exceeds the number of new equations derived, and 
this process continues as more and more equations 
are derived. This is called the closure problem 
and is the central problem of statistical turbu- 
lence theory that forces the modeling of the 
turbulence. 


justify these costs for most problems. Later in 
this paper, some examples will be presented for 
which the use of the Reynolds-stress-transport 
equations indicates some advantages that may merit 
application of these equations to certain 
problems. 

A turbulence model that takes advantage of 
the inherent characteristics of the Reynolds- 
stress-transport equations, but only adds two 
field equations to the system, is called the Alge- 
braic Stress Model (Rodi, 1976). Here the stress- 
transport equation is rewritten as 
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Levels of Turbulence Modeling 

The turbulence models required to evaluate 
the Reynolds stresses in equation (3) have been 
expressed at a variety of levels. These will be 
illustrated in this section, beginning with the 
field equations for the Reynolds stresses. By 
starting with the more complex models, it can be 
shown that some of the simpler models follow from 
the imposition of certain assumptions onto the 
more general models. 

The Reynolds-stress-transport model ( RSE) can 
be expressed as 



The terms on the right-hand side of the equation 
are divided into four quantities identified with 
the physical processes known to occur in a turbu- 
lent flow: P^j is the production of the Reynolds 

stress tensor; is the dissipation rate of the 

Reynolds stress tensor; D i j is the diffusion of 
the Reynolds stress; and 4 *^ is the pressure/ 
rate-of-strain production of J the Reynolds stress 
tensor. Of these, only the production term P i . 
can be expressed in terms of the Reynolds J 

stresses, the mean velocities, and the coordi- 
nates; consequently, it does not need to be 
modeled. The other terms are composed of third- 
order moments, moments of derivatives, or pressure 
fluctuations, all new quantities that have to be 
modeled. This level of modeling is called second- 
order closure, that is, the closure occurs in the 
Reynolds stress equations instead of in the 
Reynolds-averaged momentum equations. The equa- 
tion for e, the rate of dissipation of the turbu- 
lence kinetic energy, used to define the and 

the length scales, requires all of its termsto be 
modeled. For at least two reasons, the use of 
second-order closure is now limited to a few 
practitioners: first, the addition of seven addi- 

tional field equations increases computation costs 
by a factor of about 5 over that of the simplest 
models; second, the advantages that result from 
the use of the stress-transport equations do not 


The basic assumption of the method is represented 
by the term on the extreme right, where it is 
assumed that the convection minus the diffusion of 
the Reynolds stresses is proportional to the same 
quantities in the turbulence kinetic energy equa- 
tion. The Reynolds stresses can be written alge- 
braically as 



( 6 ) 


Recall that the production term also contains the 
Reynolds stresses, so that equation ( 6 ) has to be 
solved in an iterative manner. The kinetic energy 
and rates of dissipation are found from equations 
of the form 
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The next lower level of turbulence modeling 
utilizes the eddy-viscosity concept, that is, the 
effects of turbulence are expressed in terms of an 
effective kinematic viscosity acting on the fluid 
rate of strain. The constitutive relations 
between stress and strain that have been utilized 
are shown in the following equations from Wilcox 
and Rubes in (1980) 
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and in that from Boussinesq (1877) 
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Equation (9) accounts for the effects of the 
interaction of vorticity and the rates of strain, 
which allows experimentally observed anisotropies 
to develop in a homogeneous shear flow. Those 
effects are neglected in the classic form of 
equation (10). 

The manner of evaluating the eddy viscosity 
distinguishes most of the currently popular turbu- 
lence models. In general, the eddy viscosity can 
be written as 


v = C LV (11) 

v 

where C v is an empirically established coeffi- 
cient or function, V is a characteristic velocity 
scale, and L is a characteristic length scale. 

The various eddy-viscosity models are classi- 
fied in Table 1 ; the terminology used in the table 
will be used in the discussion of the results to 
be shown later. 

Table 1 Classification of eddy-viscosity models 


Two-equation models: 

V = yfi, L = L(k,c) 
One-equation models: 


V = yfi, L r L(y) 
One-half-equation models: 
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Zero-equation models: 


or 



L = Ky (inner zone) 


= U = C ^6* (outer zone) 
e 2 

V = V(8), L = L[(y 8)max ] 


The two-equation models utilize the kinetic 
energy and dissipation rates found from 
equation (7) and from forms of equation (8); thus, 
they require the solution of two additional field 
equations. This increases computer costs by about 
a factor of 3 over computations based on the zero- 
equation models shown at the bottom of the table. 
The reason for this increased cost is that the 
turbulence modeling equations are stiff; as a 
result, they require very small mesh dimensions 
near the surface when they are modified to allow 
their integration to the surface. It will be 
shown later that these cost penalties can be 
eliminated through the use of wall functions, that 
is, algebraic relationships that span the distance 
between the surface and the first mesh points 


located well into the regions where the turbulence 
dominates the viscous forces. 

Particular attention should be given to the 
model in Table 1 designated as a one-half-equation 
model. The "one-half" is used to emphasize the 
fact that an ordinary differential equation is 
added to the Reynolds-averaged Navier-Stokes 
system instead of another partial differential 
field equation. In the example shown here, the 
subscript 0 under the shear-stress symbol means 
an equilibrium value or one given by the zero- 
equation models indicated on the lower portion of 
the table. The subscript max means the maximum 
shear stress at some point within the boundary 
layer at the station, x. Thus, the one-half- 
equation model accounts for the lag in the devel- 
opment of the turbulence that occurs when the mean 
motion is modified by the imposition of a large 
streamwise pressure gradient. The dependent vari- 
able, the maximum shear stress in this case, is 
then used to scale the turbulence algebraically 
across the entire viscous region at station x. 

The zero-equation models (Table 1) represent 
the eddy viscosity in two zones, or layers, of the 
viscous region. Near the surface, the velocity 
s cale usually used is the local friction velocity, 
^TVp, and the length scale is merely proportional 
to the distance from the surface. In the outer 
zone, the velocity scale is the velocity at the 
edge of the viscous region and the scale depends 
on the displacement thickness (Cebeci and Smith, 
1974). In those situations where it is difficult 
to define these terms, that is, where the inviscid 
region has a non-uniform velocity field, the 
velocity and length scales have been expressed in 
terms of the local vorticity that occurs where the 
moment of vorticity is a maximum (Baldwin and 
Lomax, 1978). Since all of these quantities are 
expressible algebraically in either the dependent 
or independent variables of the mean-flow equa- 
tions, no additional partial differential equa- 
tions for field variables are required and these 
are termed zero-equation models. 


Difficulties in Applying and Developing 
Turbulence Models 

In a particular application, the level of 
turbulence model that can be considered must 
reflect the computational power available and the 
solution algorithms to be used. Limited computa- 
tional power permits only the simplest of turbu- 
lence models to be considered. Also, the turbu- 
lence model employed has a strong influence on the 
numerical behavior of the algorithms. Not all 
algorithms that can solve laminar flows can be 
applied to turbulence models, especially those 
that use field equations for turbulence quantities 
such as the Reynolds stresses, kinetic energy, or 
dissipation rate. Iteration processes are 
influenced considerably by the stiffness of the 
turbulence field equations. Also, algorithms that 
require conservative forms have difficulty with 
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the source and sink terms that are inherent in the 
turbulence field equations to reflect the turbu- 
lence production and dissipation that occur within 
the real turbulent flow. Finally, the turbulence 
models also impose constraints on the usable mesh 
configurations by requiring very close spacings 
between meshes in the vicinity of surfaces. 

This strong influence of the turbulence model 
on the numerical scheme has made it difficult to 
assess the quality of particular turbulence 
models, because often a significant part of the 
disagreement between computed results and experi- 
mental data is numerical, resulting from a pre- 
viously verified algorithm's difficulty with the 
constraints introduced by a new turbulence model. 
Once these numerical difficulties are overcome and 
the turbulence models can be applied with confi- 
dence, the resulting disagreement with experi- 
mental data for a particular flow field can be 
ascribed to the limitations of a the particular 
turbulence model being employed. At this point, 
provided fundamental experimental turbulence data 
are available to guide its improvement, it is 
possible to remodel the turbulence model to 
improve its performance for the complex flow in 
question, and to do so in a manner that does not 
alter its behavior for simpler flow fields for 
which it had already proved successful. Some 
examples of this process are shown in the follow- 
ing section where the state of the art of statis- 
tical turbulence modeling is summarized with a few 
selected flow conditions. 


APPLICATION OF VARIOUS LEVELS OF TURBULENCE 
MODELING 


In this section, in order to provide a brief 
overview of the state of the art of turbulence 
modeling, various levels of turbulence modeling 
will be applied to a variety of flow fields of 
interest to aerodynamics. Both original models 
and the consequence of "remodeling” will be shown 
in several examples. Care has been taken in most 
of these calculations to eliminate numerical 
uncertainties, so that any differences shown 
between the computed results and experimental data 
can be ascribed to turbulence-model deficiencies. 

The first flow to be considered deals with 
the effect of large streamwise pressure variations 
on the skin friction of an attached turbulent 
boundary layer (Kussoy et al . , 1978). In 
figure 2, a sketch is shown of the experimental 
arrangement. The boundary layer being measured 
was formed on the inside surface of a tube down- 
stream from a supersonic nozzle creating an axi- 
symmetric flow at M = 2.3. The pressure distri- 
bution impressed upon the boundary layer was 
created by a center body that compressed and then 
expanded the flow. The resulting pressure 
distribution, normalized by the upstream wall 
pressure where M = 2.3, is shown in the upper 
left corner of this figure. In this example, the 


pressure rise almost doubles its original value 
before relaxing. This nondimens ional pressure 
distribution remained essentially fixed as the 
wind-tunnel pressure was altered to obtain a large 
range of Reynolds numbers. The skin friction was 
measured with a buried wire gauge (Rubesin et al., 

1975) . The skin-friction distribution was 
achieved with a fixed skin-friction gauge by 
moving the centerbody upstream and downstream 
about 20 cm along its axis. This provides about 
the same distribution of skin friction that would 
be achieved with a fixed centerbody and a moving 
skin-friction gauge, because the boundary layer 
ahead of the pressure rise was close to uniform 
thickness . 

Measured and calculated skin-friction distri- 
butions are shown in figure 2 for four values of 
Reynolds number. Five levels of turbulence model- 
ing are shown: a zero-equation model by Cebeci 

and Smith (1974) in its original form and as sub- 
sequently "remodeled" to account better for 
streamwise pressure gradients (Cebeci and Smith, 
1974); a two-equation model by Wilcox and Rubesin 
(1980); and two stress-transport models. One of 
the stress-transport models utilizes a length 
scale that is prescribed algebraically (Sullivan, 

1976) , whereas the other uses a field equation for 
specific dissipation to create a length scale 
(Wilcox and Rubesin, 1980). 

With the exception of the original zero- 
equation model at the lowest Reynolds number, 
which tends to separate prematurely, all the 
models show similar variations relative to the 
experimental data. The computed results drift 
through the data with increases in the Reynolds 
number. Computed results that are low at the 
lower Reynolds numbers tend to be high when the 
Reynolds numbers are high. No one model is 
clearly superior to the others over the entire 
Reynolds-number range, though they do seem to 
maintain their relative positions to each other. 
One can conclude from this work that for an 
attached boundary layer there is little advantage 
to using the more complex models, but if it were 
necessary to know the skin friction in attached 
boundary layers in strong streamwise pressure 
gradients to a few percent, improvements to all 
the levels of models would be required. 

Streamwise curvature of a surface has a pro- 
nounced effect on the behavior of a boundary layer 
flowing over it (Bradshaw, 1973). This is demon- 
strated in figure 3, where skin friction and 
shape-factor data obtained on a convex surface in 
an adverse pressure are shown (So and Mellor, 

1972). The data are expressed in terms of the 
streamwise distribution of the skin-friction coef- 
ficient and the shape factor of the turbulent 
boundary layer flow. Also shown in the figure are 
the computed results utilizing four turbulence 
models (Wilcox and Rubesin, 1980). The dashed and 
dot-dashed curves are based on a two-equation and 
a Reynolds-stress transport model in which the 
effects of curvature have been neglected. These 
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skin-friction coefficient curves depart consider- 
ably from the experimental data. The shape factor 
on a flat plate would remain at a value of about 
1.4, again showing curvature effects that are 
large in affecting the boundary-layer velocity 
profiles. The solid and dotted curves represent 
the same models modified to account for the sur- 
face curvature through the introduction of a 
curved coordinate system. The modification to the 
two-equation model also required the introduction 
of a new centrifugal force term in the turbulent 
kinetic energy equation; that term vanishes on a 
flat surface. This was done in a rather ad hoc 
manner, utilizing some guidance from the the form 
of the Reynolds-stress equation for v'v’. On the 
other hand, no physical modifications were made in 
the Reynolds-stress model. In both of these 
models, the specific rate-of-dissipation equation 
used to define the length scale was not altered 
except for the coordinate modifications. Figure 3 
shows that the effects of the streamwise curvature 
in these data are represented very well by these 
latter computations. 

It is very significant that no modifications 
had to be made to the Reynolds-stress turbulence 
model to allow it to capture the effects of 
streamwise curvature; it is thereby an example in 
which, at least in principle, the more complex 
form of modeling has a fundamental advantage. 
Another conclusion, however, is that the simpler 
two-equation model can also be "remodeled" to give 
good results on a convex surface, without altering 
its behavior, for cases in which the streamwise 
curvature is absent. A cautionary note is in 
order here: concave surfaces may not be calcu- 

lated well with any of these models because the 
models completely ignore the possibility of the 
existence of Gortler vortices that are present 
over such surfaces and that can enhance the trans- 
port within the boundary layer. 

Figure 4 shows the experimental data and 
computed results for the distribution of surface 
pressure and skin friction in the vicinity of a 
normal shock wave in an M = 1.48 airstream. 

These measurements were obtained on the inner 
surface of a circular tube with a static pressure 
tap and a flush hot-wire skin friction gage 
(Mateer et al., 1976). The shock-wave position 
could be moved relative to the fixed gauges by 
increasing or decreasing the blockage of a down- 
stream shock generator. 

The experimental data are compared with com- 
putations utilizing four different turbulence 
models: 1) two-equation models from Wilcox and 

Rubesin (1980), Jones and Launder (1972), and 
Chien (1982); and 2) a zero-equation model from 
Cebeci and Smith (1974), which did not incorporate 
the pressure-gradient modification discussed in 
connection with figure 2. All of these models 
were used in computations that extended to the 
surface. It should be noted that each of the 
models performs quite well on a flat plate at 
uniform pressure. The models of Jones and Launder 


and of Chien are very similar, except that the 
latter was more closely fitted to experimental 
channel data. 

When these models are applied to the shock- 
wave and boundary-layer interaction of this exper- 
iment, it is found that all of them yield excel- 
lent surface-pressure distributions. The skin- 
friction results, however, are far from satisfac- 
tory for all but the Wilcox-Rubesin model. 

The Chien model, which behaves so well on a flat 
plate, yields surprisingly poor results in the 
vicinity of the shock wave. The Jones-Launder and 
Cebeci-Smith models fail to agree, in opposite 
ways, with the skin friction, the latter model 
indicating an extensive region of separation where 
none existed experimentally. 

Qualitatively, this behavior is consistent 
with the behavior of the Cebeci-Smith model in 
figure 2, when it was not altered to account for 
pressure gradients. From a comparison of the 
different skin-friction results based on the dif- 
ferent two-equation models, it can be concluded 
that models at the same level that are "cali- 
brated" with flat-plate data can behave much dif- 
ferently when applied to a complex flow. This 
suggests that in general, the predictability of 
any turbulence model should be suspect until it is 
verified for a variety of complex flows. 

In an effort to determine why the Jones and 
Launder and the Chien turbulence models predicted 
such different skin-friction results, Viegas and 
Rubesin (1983) eliminated the near-surface model- 
ing in these turbulence models by applying the 
concept of wall functions to the computations 
shown in figure 4. This concept requires the 
first mesh point off the surface to lie well 
within the fully turbulent flow and therefore 
connect the flow at this point to the surface 
conditions with algebraic formulas. For example, 
a wall function in its simplest form is the con- 
ventional "law of the wall" applied to a zero- 
equation turbulence model on a flat plate with 
small streamwise pressure gradients. The two- 
equation models require more complex wall func- 
tions that account for the growth of turbulence 
kinetic energy in regions of strong pressure gra- 
dients, the mean values of kinetic energy and 
dissipation rates within the first mesh volume off 
the surface, and the changes to the law of the 
wall caused by local pressure gradients. 

Figure 5 shows the results of applying such 
wall functions to the test conditions described in 
figure 4. Again, the surface pressure is well 
predicted when the different models use wall func- 
tions at distances from the wall shown by the 
corresponding upstream values of y+ indicated in 
the figure. The skin friction results show 
clearly that it is the individual near-wall treat- 
ments of the two-equation models by Jones and 
Launder and by Chien that causes their anomalous 
behavior shown in figure 4. 
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With wall functions, the results given by the 
two models are essentially, identical and quite 
consistent with the other two-equation model bv 
Wilcox and Rubesin. The Cebeci -Smith model with 
its wall function is quite consistent with the 
earlier results where' it was computed to the sur- 
face, except in the immediate region of the sepa- 
ration zone where the pressure gradients are suf- 
ficiently high to bring into question the use of a 
usual law of the wall. It is gratifying that the 
wall functions not only collapse the results of 
the different two-equation models, but that 
results agree quite well with the experimental 
skin-friction data. Another significant bonus of 
the use of wall functions is that the costs of the 
computations were reduced by a' factor of 8 because 
of the reduced number of mesh points needed to 
resolve the near-wall region and the elimination 
of much of the stiffness of the turbulence model- 
ing equations: In this case, accuracy and economy 

were achieved together. 

Figure 6 shows the results of improvements to 
two different levels of turbulence models applied 
to the computation of a shock-wave and boundary- 
layer interaction in transonic flow. The experi- 
mental data used in this comparison are from 
Johnson et al. (1982) and they were measured on an 
ax isymmetric circular-arc bump Indicated schemati- 
cally in figure 6. The free-sftream Mach number of 
the examples shown here is M = 0.876 and the 
unit Reynolds number is 13-1 * 106°/m. 

The upper part of figure 6 shows the distri- 
butions of surface pressure measured in the exper- 
iment and computed from four turbulent models. 

The abscissa begins at mid-chord of the bump and 
continues about a half-chord downstream after its 
trailing edge. The lower part of the figure shows 
the distribution of skin friction from two of the 
computations and the experimental points' of sepa- 
ration and attachment, as they were indicated by 
an oil-film technique. This demonstrates that 
there was a rather extensive region of separated 
flow in this experiment. 

The computed surface pressures based on the 
zero-equation turbulence model of Cebeci and Smith 
(1974) are represented by the dashed line. These 
results indicate that this model predicts a shock- 
wave position that, is well downstream of its mea- 
sured location. The computations also show much 
higher pressure in the region of separation than 
do the data. The dotted line represents the com- 
puted results obtained with the half-equation 
model of Johnson and King (1985). The latter 
model accounts for the history of the development 
of maximum shear stress in the streamwise direc- 
tion through the use of an ordinary differential 
equation. This maximum shear stress is then used 
to scale an algebraic model , similar to that of 
Cebeci and Smith, at each station. Although these 
experimental data, for a range of Mach and 
Reynolds numbers, were used to evaluate the lag 
constants used in this model, the resulting agree- 
ment of the computed results with the shape of the 


entire pressure distribution for this particular 
case is remarkable. Without altering' the modeling 
coefficients’, good results have been obtained. by 
Johnson and King in comparisons with other 
two-dimensional data Under widely different flow 
conditions (Simpson etal., 1 98 1 ) . In addition, 
this model is economical to use, for it requires 
very little more computer time than do the stan- 
dard zero-equation models. 

When a two-equation model is -applied to this 
flow, the dot-dashed line (fig. 6) results for the 
predicted pressure distribution. Although it 
performs better than the zero-equation model, the 
pressure distribution calculated with the two- 
equation model also shows a shock-wave location 
that is downstream of the experimental position, 
and the pressure is overpredicfced to some extent 
in the' region of separation. Wh"en improved wail 
functions ( Viegas et al . , 1985; Rubesin and : 
Viegas, 1985) are applied to this Johes-Launder 
model, the calculated pressure distribution is 
considerably improved, giving results equivalent 
to the Johnson-King model except at the upstream 
foot of the shock wave. It should be noted that 
the fohms of the wall functions or the original 
model were not altered to "fit these particular 
pressure data. The use of wall functions with the 
Jones-Launder model also improves the prediction 
of the location of separation and reattachment as 
is seen in the skin-friction coefficient figure. 

In conclusion, figure 6 shows that two turbu- 
lence models are available that can give good 
surface-pressure distribution results in a tran- 
sonic flow over an airfoil-like body. For two- 
dimensional flows with a closed separation region, 
the model of choice is the Johnson-King model ^for 
accuracy and economy of application. It is not 
clear at present, however, how or if this model 
can be extended to three-dimensional flows. The 
two-equation model is directly extendable to three 
dimensions; however, it requires some additional 
modeling regarding the degree of isotropy of the 
eddy viscosity, and the wall-function approach in 
its present form is constrained to no skewing of 
the mean-velocity vector within the first mesh 
volume’ off the Surface. 

Figure 7 shows the results of applying sev- 
eral levels df turbulence modeling (Coakley, 1986) 
to the prediction of the lift and drag of an 
RAE 2822 airfoil for three test Mach numbers and 
angles of attack (Cook et al., 1979). The com- 
putations are shown here in lift-drag form to 
allow comparison of the results relative to a 
constant lift-drag ratio as shown by the lines 
passing through the three starred data points. 
Coakley (1986) used six different turbulence 
models: two zero-equation models, that of Cebeci 

and Smith and of Baldwin and Lomax; the Johnson- 
King half-equation model; and three versions of 
the two-equation model, Chien's and two originated 
by Coakley. The symbols used for the results 
found from these models are indicated on the 
figure . 
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The comparison of the experimental data and 
the computed results can be observed from two 
points of view. First, a small vector distance 
between the experimental and computed points 
reflects good agreement between the experimental 
data and computation for both 1) surface-pressure 
distributions and shock-wave locations, and 
2) skin-friction distributions. Second, the 
alignment with the constant L/D line suggests that 
a particular model can yield a result that happens 
to provide good L/D, although in reality it misses 
the details. 

The cases represented by the three Mach num- 
bers correspond, with increasing Mach number, to 
unseparated subcritical flow, unseparated super- 
critical flow, and separated supercritical flow, 
respectively. Except for two cases at M = 0.73, 
the computed results yield values of lift and drag 
that are higher than the measurements. At all 
three Mach numbers, the Johnson-King model lies 
closest to the data, even at the severest flow 
with separation. The other models generally 
become poorer as the flow becomes more complex. 

The q-w2 model of Coakley, while second-best 
overall in agreeing with the experimental results, 
yields a much poorer lift-drag slope than do the 
other models when separation takes place. 

The earlier conclusion that the Johnson-King 
model is the model of choice for two-dimensional 
transonic flows is further supported by these 
data. For those interested only in the ratio of 
lift to drag, the simplest models of Cebeci-Smith 
and Baldwin-Lomax are reasonably adequate. These 
conclusions, however, are subject to the caveat 
that these data resulted from flows that were 
entirely attached or flows that separated and then 
reattached on the airfoil. Cases with massive 
separation may favor other turbulence models. 

The next few figures will show the perfor- 
mance of various turbulence models for three- 
dimensional flows. Although the first two 
examples are strictly two-dimensional from a 
mathematical viewpoint, that is, axisymmetric flow 
or flow over an infinitely long, swept-back air- 
foil, the flows exhibit the skewing of the veloc- 
ity vector parallel to the surface that is char- 
acteristic of three-dimensional boundary-layer 
flows . 

Figure 8 describes skin-friction results 
obtained on the surface of a circular cylinder 
whose axis was mounted col inear with that of a 
low-speed wind tunnel (Higuchi and Rubesin, 1 98 1 ; 
Driver and Hebbar, 1985). The particular feature 
of this experiment was that the cylinder was seg- 
mented so that a central portion of the cylinder 
could be rotated to introduce a shear-driven cross 
flow (see sketch in fig. 8). In the stationary 
region ahead of the rotating segment, the boundary 
layer deve lops the four Reynolds stresses, 
u'u', v'v', w'w', and u'v', which are characteris- 
tic of two-dimensional flow. Over the spinning 
portion, the additional Reynolds stresses, 


u'w' and v’w', are activated, and it is the latter 
stress that causes the mean-velocity vector to 
skew in a three-dimensional boundary layer. In 
the experiment, the relaxation of the Reynolds 
stresses back to their two-dimensional character 
is studied on the downstream stationary cylinder. 

Figure 8 shows the development of the axial 
and cross-flow skin friction with distance along 
the stationary downstream cylinder. The experi- 
mental data were measured with a variety of tech- 
niques: buried hot-wire gauges, surface fences 

(Higuchi, 1983), and momentum integral balances. 
Computed results based on three levels of turbu- 
lence modeling are also shown. These include the 
zero-equation Aquilar (1976) model which is an 
extension to three dimensions of the Cebeci-Smith 
model, and the two-equation and Reynolds-stress- 
transport models of Wilcox and Rubesin. Note that 
both the zero-equation and the two-equation models 
utilize the additional assumption that the eddy 
viscosity is a scalar quantity that acts equally 
on the axial and cross-flow strains. Generally, 
all the models capture the behavior of the relaxa- 
tion of the skin friction in both directions. Of 
the three models, however, the Reynolds-stress- 
transport model yields the best agreement with the 
data, except for the short region at the beginning 
of the stationary cylinder where the relaxation 
process is occurring very rapidly. Since this 
flow remains attached, it is believed the short- 
comings of the zero-equation and two-equation 
models, both of which are eddy-viscosity models, 
are caused by the the scalar eddy-viscosity 
assumptions and not by problems of length scale. 


Figure 9 shows oil-film streaks on the sur- 
face of flow over a transonic wing, swept back at 
an angle of 32° in an airstream at M = 0.74 
and Rec = 4.7 * 10^ (Mateer and Brosh, 1983). 

The airfoil shape is a NACA 0012 at zero angle of 
attack. The wing spans the wind tunnel, but the 
experiment is unique in that the side walls of the 
wind tunnel were contoured to follow the expected 
inviscid flow, thereby causing the model to behave 
as one of infinite length where surface pressures 
are essentially uniform in the spanwise direction. 
Also shown on the figure are the surface stream- 
lines as computed by three levels of models. 

Line 1 refers to a Reynolds-stress-transport 
model; line 2 to a two-equation model; and lines 3 
and 4 refer to a mixing-length model. These 
models are identical to those shown in the pre- 
vious figure for the spinning cylinder experiment. 


Two results for the mixing-length model are 
shown to demonstrate the influence of the location 
of transition, which can be specified in this 
turbulence model. In the other models the loca- 
tion of transition is established by stability 
criteria built into the models. The mixing-length 
model shows an increased deviation. from the oil- 
flow data as transition is moved downstream. As 
with the spinning cylinder, the Reynolds-stress- 
transport model gives slightly better results as 


229 



represented by the oil-flow pattern near the suc- 
tion peak. 

Figure 10 also shows a comparison of experi- 
mental and computed surface skin-friction lines 
for a three-dimensional flow (Horstman et al., 
1985). In this case, the test body consists of a 
cylinder followed by skewed flare. The airstream 
conditions are M = 3 and the length Reynolds 
number of the boundary layer just upstream of the 
flare is 18 x 10^. The flare is a 60° total-angle 
cone tilted to an angle of attack of 23°. The 
computations for this flow were based on the 
Jones-Launder two-equation turbulence model with 
boundary conditions supplied by wall functions 
(Viegas and Rubesin, 1983). The use of wall func- 
tions was critical to reduce the cost of these 
computations . 

The oil-flow lines on the left of this figure 
show that for these conditions the shock wave 
generated by the flare is sufficiently strong to 
cause the boundary layer on the cylinder to sepa- 
rate well upstream of the cylinder-flare junction 
and to reattach about halfway up the flare. Con- 
siderable amounts of turning of the near-surface 
flow takes place within this separation zone and 
in the immediate region after reattachment. The 
computed surface shear lines show the general 
character of measured oil streaks; however, the 
details are off to quite an extent. The computed 
upstream separation location is only about one- 
half its experimental distance from the flare 
junction. The poor behavior of the computations 
is also reflected in a generally smaller upstream 
movement in the computed results in the vicinity 
of the flare. Although the reattachment region is 
predicted fairly well, the flow direction on the 
cylinder downstream of the flare also departs from 
the experimental data. 

How these differences are felt quantitatively 
is shown in a comparison of the computations with 
the experimentally measured surface pressures in 
figure 11, where the surface pressures on the 
windward ray are plotted along the cylinder, 
flare, and afterbody. The smaller predicted sepa- 
ration zone is again evident in this figure. 

These quantitative results of the computations 
show the position of the maximum pressure to be a 
bit upstream of the measured location and to yield 
a smaller value than measured. These results 
suggest that the computed side force and pressure 
drag are low. Recently, Brown et al . ( 1 987 ) made 
measurements of mean flow and the Reynolds 
stresses under these same flow conditions. These 
data have not yet been analyzed with the view of 
improving the turbulence modeling; however, this 
activity is under way and should lead to a 
remodeling of the turbulence model that will 
result in computations that are considerably 
improved over those shown here. 


CONCLUDING REMARKS 


Certain observations can be drawn from this 
brief survey of recent developments in statistical 
turbulence modeling applicable to aerodynamics. 

For attached boundary layers upstream of pressure 
gradients, all levels of turbulence modeling give 
about the same results. This is not surprising 
because the attached flat-plate boundary layer has 
been used universally to "calibrate" these models. 
Under severe pressure gradients, however, many 
popular models still yield results that differ 
from experimental skin-friction data by as much as 
k0% when the Reynolds numbers are about 100 x 10° 
or larger. This is a turbulence modeling topic 
that requires additional attention. 

The manner of modeling streamwise curvature 
effects is also still an open question. The sim- 
plest zero-equation models use an ad hoc correc- 
tion method recommended by Bradshaw (1973). On a 
convex surface, the two-equation models account 
for curvature by either increasing the rate of 
production of dissipation (Ha and Lakshminarayana, 
1980), by diminishing the rate of destruction of 
dissipation (Launder et al., 1977), or by letting 
the dissipation rate alone and altering the rate 
of production of kinetic energy by introducing the 
work of centrifugal forces (Wilcox and Rubesin, 
1980 ; Pulliam et al. 1985). The full Reynolds- 
transport model, on the other hand, requires no 
changes other than the introduction of curved 
coordinates (Wilcox and Rubesin, 1980). This 
variety of approaches for the introduction of the 
effects of surface curvature into the zero- and 
two-equation turbulence models illustrates the 
fact that these models are merely such gross 
approximations to the physical turbulence that 
there is no unique way of modeling a particular 
phenomenon. Consequently, a model can only be 
gauged by its success or failure in its applica- 
tion to a particular flow field or to a range of 
flow fields. There is no a priori way of asses- 
sing the accuracy and breadth of application of 
turbulence models other than by comparison with 
experimental data. From such comparisons, it is 
generally observed that the higher levels of tur- 
bulence models tend to have broader ranges of 
application, but that for a given application they 
are likely to be less accurate than simpler models 
that have been fitted to the particular applica- 
tion. The simple, or low-level, models, tnough 
computationally economical, need a data base for 
every class of flow likely to be encountered, and 
their use places an enormous burden on the numbers 
of experiments required. Without such a data 
base, the simple models are not much sounder than 
intelligent guesses. It is this burden on experi- 
ment that stimulates the search for more general 
models that are sufficiently accurate for a vari- 
ety of applications even though they may cost more 
to compute. The user must gauge these trade-offs 
in deciding which models he will employ. 
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With respect to the rate of progress in tur- 
bulence modeling, the development of the Johnson- 
King model is a good example. From the time of 
the first experiment on airfoils in which Johnson 
was measuring turbulence quantities (Johnson and 
Bachalo, 1980) to the recent successful tests of 
the model he and a colleague developed (Johnson 
and King, 1985), 6 years elapsed. Much of this 
time was spent by Johnson in analyzing his data 
with a variety of borrowed computer codes, some of 
which were in development and were extremely dif- 
ficult and time-consuming to operate. During this 
period, he also performed an additional experiment 
on the axisymmetric bump model that was used to 
complement the original airfoil data. This short 
history demonstrates that turbulence modeling is 
composed of a variety of activities that encompass 
analysis of the modeling equations, experimenta- 
tion, the use of and modification of complex com- 
puter codes, and, most important, a dedication to 
do the necessary labor to improve an existing 
model . 

When one or two individuals attempt this, the 
level of effort required is bound to take consid- 
erable elapsed time. Something as relatively 
direct as developing wall functions for compres- 
sible flows, utilizing the experimental work of 
others and making small modifications to existing 
computer codes (Viegas and Rubesin, 1983; Viegas 
et al., 1985), took a little longer than 2 
years. Because of its dependence on a variety of 
disciplines, turbulence modeling is a slow process 
and its development can be accelerated only by an 
increased level of the coordinated efforts of 
dedicated people possessing a variety of 
complementary talents in analysis, 
experimentation, and computational-code develop- 
ment. This is especially required now as the CFD 
community begins to emphasize complex three- 
dimensional flow fields and requires improved 
turbulence models for such flows. 
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Fig. 2 Compressible boundary layer skin friction 
in an axial pressure gradient. 
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Fig. 3 Boundary layer development on a convex Fig. 4 Distribution of surface properties in the 

surface in a low speed flow with an axial pressure vicinity of normal shock wave-computations with 

gradient, (a) Skin friction; (b) shape factor. integration to the wail. 


235 





(x - x 0 )/S 0 


Fig. 5 Distribution of surface properties in the 
vicinity of normal shock wave-computations with 
wall functions. 
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Fig. 6 Distribution of surface properties on an 
axisymmetric bump in transonic flow. 


236 




ORIGINAL PAGE IS 
OF POOR QUALITY 


-(> EXP. DATA 
o CEBECI-SMITH (0-EQ) 

□ BALDWIN-LOMAX (0-EQ) 
O JOHNSON-KING O/ 2 -EQ) 
LINES OF CONST L/D 


2-EQ MODELS 
▲ CHIEN (k-e) 

A COAKLEV (q-cjl) 
M COAKLEV (q-w2) 




Fig. 7 Lift-drag relationships for an RAE 2822 
airfoil . 


Fig. 8 Skin friction in the relaxation zone 
downstream of transversely sheared boundary layer. 



Fig. 9 Comparison of experimental and computed 
surface flow lines on a swept, infinitely long 
wing-NASA 0012 at zero angle of attack. 
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